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We solve the 5=1/2 infinite-range random Heisenberg Hamiltonian in the paramagnetic phase 
using quantum Monte Carlo and analytical techniques. We find that the spin-glass susceptibility 
diverges at a finite temperature T g which demonstrates the existence of a low-temperature ordered 
phase. Quantum fluctuations reduce the critical temperature and the effective Curie constant with 
respect to their classical values. They also give rise to a redistribution of spectral weight in the 
dynamic structure factor in the paramagnetic phase. As the temperature decreases the spectrum 
of magnetic excitations gradually splits into quasi-elastic and inelastic contributions whose weights 
scale as S 2 and S at low temperature. 
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It is well known that the combined effects of ran- 
domness and frustration may lead to spin glass behav- 
ior in classical disordered magnets at low temperature. 
The most widely studied spin-glass Hamiltonian is the 
Sherrington-Kirkpatrick model However, only a 

small fraction of the vast amount of work Q devoted to 
this model directly addresses the role of quantum fluc- 
tuations. In a notable early paper Bray and Moore ||] 
first formulated the theory of the quantum Sherrington- 
Kirkpatrick model and showed that it reduces exactly to 
an effective single site problem in imaginary time. Using 
a variational approach, these authors argued that a spin 
glass ordered phase occurs below a finite critical temper- 
ature for all values of S. Much more recently, Sachdev 
and Ye Q discussed a generalized spin-glass Hamiltonian 
in which the spin components become the generators of 
the group SU(M) and the states span a representation of 
the group labeled by an integer rib {M = 2 and rib = 2S 
for physical spin-5 operators). The model can be solved 
exactly when M and — > oo with k = n^/M finite. In 
this limit, Sachdev and Ye argued that the ground state 
is either a spin glass for large values of K (that plays the 
role of S), or a spin fluid below a critical value n c . In the 
spin fluid phase the local dynamic susceptibility exhibits 
unconventional behavior ~ ln(— ) at T = Q. 

In view of these results, the nature of the ground state 
of the quantum Sherrington-Kirkpatrick model remained 
an open problem, the main question being whether quan- 
tum fluctuations for low enough S may prevent the in- 
stability towards spin-glass order present in the classical 
case ||. 

In this paper we answer this question by means of an 
exact numerical solution of the 5=1/2 model using a 
quantum Monte Carlo technique, and obtain analytical 
expressions for the asymptotic forms of the imaginary 
part of the local dynamic spin susceptibility x"( w ) m the 
limits T — > and T — > oo. We find that the paramag- 
netic solution is unstable towards spin glass order at a 
finite temperature T g . The transition temperature and 



the effective Curie constant are reduced with respect to 
their classical values by quantum effects. The analysis 
of the paramagnetic dynamic correlation function shows 
that the spectrum of magnetic excitations splits at low 
temperatures into two well defined contributions. The 
first one is similar to the low-frequency dynamic response 
function of classical fluctuating paramagnets ||. The 
second one is an inelastic contribution at higher frequen- 
cies, characterized by a large energy scale. At high T, 
all the spectral weight S(S + 1) is concentrated in the 
quasi-elastic feature. With decreasing temperature part 
of the weight is transferred from low to high energies un- 
til, when T — ► 0, the intensities of the quasi-elastic and 
inelastic parts reach the asymptotic values S 2 and S, re- 
spectively. 

The Sherrington-Kirkpatrick Hamiltonian is, 
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where Si is a three dimensional spin-1/2 operator at i— th 
site of a lattice of size N. The exchange interactions 
Jij are independent random variables with a gaussian 



distribution with zero mean and variance J = (Jfj) 1 ^ 2 - 

In the limit TV — > oo, the free-energy per spin T has 
been derived by Bray and Moore || using the replica 
method to do the average over the quenched disorder 
and a Hubbard-Stratonovich transformation to decouple 
the different sites, 
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Here, T is the time-ordering operator along the 
imaginary-time axis, < r < /?, and the trace is taken 
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over the eigenstates of the spin-1/2 operator. The local 
dynamic correlation function in imaginary-time Q(t), is 
determined by functional minimization of (0), 
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Q(r) = -<T5(r) • 5(0)), 
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main source of error in the results comes from the statis- 
tical noise due to the random sampling. The efficiency of 
our code allowed us to sensibly reduce it by performing 
10 5 full updates in the last iteration. It is remarkable 
that we faced no "sign problem" even down to the lowest 
temperature considered, T = 0.02 J. 



where the thermal average is taken with respect to the 
probability associated to the local partition function, 



Zioc = TrT exp 
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The spin-glass transition temperature may be calcu- 
lated from local quantities. It results from the instability 
criterion H Jxioc = 1, with the static local susceptibility 

Xioc = Jq drQ{r). 

In order to set up a numerical scheme for the solution of 
the model we find it necessary to perform an additional 
Hubbard-Stratonovich transformation and rewrite Zi oc 
as 



Zioc = J Vf} exp 
x TrT exp 
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Zioc is the average partition function of a spin in an ef- 
fective "time" -dependent random magnetic field Jff(r) 
distributed with a gaussian probability. This formula- 
tion is well suited to implement a quantum Monte Carlo 
algorithm. The imaginary time axis is discretized into L 
time slices and the time-ordered exponential under the 
trace in Eq. ^| is written as the product of L matrices of 
2x2 using Trotter's formula. We performed calculations 
for ft < 50 and L < 128 (keeping JAr = J/3/L < 0.5). 
There are two important technical remarks: firstly, it is 
crucial that each trajectory ff(r) and its time-reversed 
partner rf(/3 — r) be considered simultaneously in order 
to obtain a real probability measure M. Secondly, the 
sorting procedure is formulated in the frequency domain: 
the phase space for the simulation consist of all the real- 
izations of fj(oj n ) with the integer n labelling the bosonic 
Matsubara frequencies. This new set of variables presents 
the advantage of being much less correlated than the 
original one. An elemental Monte Carlo move is thus 
to propose a change in the complex field rji(u! n ) for given 
i = x,y,z and n; a full update of the system is com- 
pleted after elemental moves have been attempted for all 
directions and frequencies. The numerical procedure is 
as follows: i) an initial Q{ui n ) is used as input in (|^). 
ii) the spin-spin correlation function is obtained using 
Monte Carlo, iii) a new Q(ui n ) is calculated from the 
self-consistency condition (^) and used as a new input in 
step i). This procedure is repeated until convergence is 
attained which typically occurs after 5 iterations. The 
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FIG. 1. Q(t) as function of t//3 for several values of /3J . 
The crosses correspond to QMC data. The error bars are well 
below the size of the crosses. Solid lines are fits obtained with 
the function Q pa r[T) described in the text. 

In Fig. |l] we show the correlation function obtained 
with this method which exhibits the following qualita- 
tive features. For T >> J, Q(t) is nearly constant and 
close to S(S + l)/3 = 1/4, its classical limit. In con- 
trast, for T « J, Q{t) rapidly decreases at both ends 
of the interval [0, 0\ and then varies slowly remaining 
near the value S 2 /3. The behavior in the intermediate 
temperature range is a smooth interpolation between the 
two extreme cases. We obtain xioc by numerically inte- 
grating Q{t) and find that it crosses over between two 
limiting forms, \ioc ~ S(S + 1)/(3T) for T » J and 
Xioc ~ S 2 /(3T) for T « J, implied by the asymptotic 
behavior of Q(t) described above. The susceptibility 
thus obeys a Curie law at both ends of the temperature 
range with the low T effective Curie constant reduced 
by a factor of 3 with respect to its classical value. Note 
that the low- temperature behavior of xioc{T) differs con- 
siderably from the prediction of the large-M model Q] 
which is Xioc(T) ~ ln(l/T) for small k. On the other 
hand, the full temperature dependence of xioc is remark- 
ably close to that predicted by the variational approach 
This is easily understood at high and low temper- 
atures where the variational ansatz of Bray and Moore 
@! Q( T ) = const, represents well the actual numerical 
solution (cf. FigJ^). However, the close agreement at in- 
termediate temperatures is rather unexpected. The cri- 
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terion for the appearance of spin glass order, Jxioc = 1, 
is fulfilled at T g w 0.142 J that compared to the classi- 
cal value T g — 0.25 J shows the importance of quantum 
fluctuations in this system. 

We now turn to the discussion of the frequency depen- 
dence of the local susceptibility. This involves extracting 
real-frequency information from results obtained on the 
imaginary axis, which is in general a rather difficult task. 
However, it will be seen below that in the present case 
the dissipative part of the local response Xioci 1 ^) mav be 
found analytically at low and high temperature. This al- 
lows us to construct an interpolating function that accu- 
rately describes our numerical data throughout the whole 
temperature range, enabling us to extract physical infor- 
mation from the imaginary-time Monte Carlo results. 

The dynamics of the system for T >> J is on general 
grounds expected to be controlled by a single relaxation 
rate ujl = O(J). This assumption implies 
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where the relaxation function F is constrained to be nor- 
malized to one and have a finite second moment || . The 
simplest function fulfilling these conditions is a gaussian. 
With this ansatz, the response is completely determined 
by the sum-rule 
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coxLM = 2J 2 / Q 2 (r), 
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that is derived using the generic /-sum rule for spin 
systems and It then follows from (^) and (Q) that 
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with = 2J 2 S(S+l)/3[l-S(S + l)(f3J) 2 /18} forT » 
J. It can be shown (?]] that Eq. || reproduces the first few 
orders of the high-temperature expansion of Q(t). This 
expression is thus correct for T >> J. 

With decreasing temperature, the assumption of a sin- 
gle relaxation rate breaks down. At T < < J the existence 
of two different characteristic times suggested by the nu- 
merical data of Fig. |l| must be reflected in the emer- 
gence of well separated energy scales in the frequency 
dependent dynamic response. Indeed, this behavior fol- 
lows from an approximate analytical solution of Eqs. ^ 
and H that can be shown to be exact in the T — » limit. 
From Eq. [| the problem can be thought of as that of a 
single spin in a fluctuating effective magnetic field Jrj(r) . 
We will see below that at low T this effective field is 
dominated by its uj n — component. Therefore, it is 
convenient to split the effective field into a constant part 
ho = Jff(u> n = 0)/t/p and a small r-dependent part with 
\dh(T)\ « \ho\. We thus write (g) 



Q(r) = r (Ql(t, Mk + 2(Qt(t, fc )>K 
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where Ql and Qt are, respectively, the longitudinal and 
transverse response functions in an applied field ho hav- 
ing formally integrated out the fields 6h(r). The angular 
brackets denote the average with respect to the isotropic 
distribution P(|/ig|). For a given ho, the imaginary part 
of the transverse response function x'ri 10 ) has a peak at 
u) = \ho\ whose width T is proportional to the square of 
the amplitude of the fluctuating field at the resonance 
frequency. A simple estimate || shows that for T « J, 
P(|/io|) is maximum at \ho\ = wt ~ J 2 Sxioc, which is 
large at low temperature. Assuming for the moment that 
|5/i(r)| << \ho\ (i.e., setting Y = 0) one can perform the 
average over P(|/io|) and obtain an approximation for the 
dissipative part of (Qt(t, ho))^ Q , 
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as 



Using this equation and the fluctuation-dissipation 
theorem we estimate Y/ujt = 0(T/J), showing that the 
assumption leading to ([lO]) is indeed correct. Therefore, 
the expression above is asymptotically exact for T — > 0. 

The high-frequency scale ojt ~ J 2 /T where Xt( w ) 
is sharply peaked, is associated to the initial rapid de- 
crease in Q(t) observed in our low-T simulations. The 
remaining contribution to the dynamic response func- 
tion x'Li 10 ) comes from the relaxation of the longitudi- 
nal magnetization. As the amplitude of the fluctuating 
field |£/i(t)| << \ho\ we expect this process to be slow. 
Its frequency ujl may be found irrespectively of the de- 
tailed form of Xz,( w ) using Eqs. 0, || and [l(] which yield 
ujl oc T 2 0| . This small energy scale is associated to 
the slowly varying part of Q(r) observed at low temper- 
ature. Using again the ansatz (g) for y?(w), the sum-rule 
(0) leads to an expression similar to (|8[), except that the 
prefactor of the exponential is now simply (3S 2 /'i. Notice 
that the dynamics of magnetic fluctuations that emerges 
from the above arguments bears no resemblance to that 
of the spin-fluid state of the large-M model. 

Comparison between the low- and high-temperature 
results implies that when T decreases, a fraction of the 
spectral weight of the quasi-elastic peak at small-cj is 
transferred to the high energy excitations described by 
Eq. [l^. This redistribution of intensity, closely related 
to the reduction of the effective Curie constant discussed 
above, is a distinctive quantum effect: the strength of 
the inelastic feature relative to that of the quasi-elastic 
peak is 0(1/ S) and vanishes in the large-S 1 limit. The 
analytical results just discussed suggest a parametriza- 
tion of the spin-spin correlation function Q[t) that con- 
tains the exact asymptotic forms at both high and low 
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T and smoothly interpolates between them. The pro- 
posed interpolation function is defined as Q P ar( T ) = 
S(S+1)/3\p$ l (t) + (1-p)$t{t)] where, 
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where and $t are the imaginary-time equivalents of 
Eqs. ^ and |l0| normalized such that < &l,t(0) = 1, and u>l 
and lot are now parameters corresponding to the width of 
the central peak and the characteristic scale of the high- 
energy excitations, respectively. The third parameter, p, 
controls the transfer of spectral weight between the two 
components of the magnetic response. At high T there 
is a single energy scale and only quasi-elastic intensity 
is present as p — > 1, while at low T inelastic intensity 
appears and p — > S/(S + 1) its lower bound. 

Using this expression we obtain highly accurate fits of 
our numerical results as demonstrated in Fig. [l]. While 
Qpar(T~) contains the correct high and low T limits by 
construction, it is remarkable that the quality of the 
agreement remains excellent through all the temperature 
range. This gives us confidence that the essential physics 
is indeed captured by our parametrization of Q(t). It 
is now interesting to go back to real frequency and plot 
x" oc (o;)/cj which is an experimentally accessible quan- 
tity. The results shown in Fig. ^ illustrate the gradual 
emergence of the high-frequency features in the dynamic 
response. 
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FIG. 2. The relaxation function F(ui) — \ '(ui) / '(nxiocu) 
for j3J = 3 (long dashed), 7 (dashed) and 14 (solid). Inset: 
T— dependence of the parameters p (circles) and u>l (squares) 
that have the limiting behavior predicted by the theory. 

Our results are strictly valid above T g , since below this 
temperature the paramagnetic state is unstable . Nev- 



ertheless, the study of the solution in the whole T— range 
is justified as it should be kept in mind that the precise 
value of T g depends on the details of the model Hamil- 
tonian. For instance, lowering the dimensionality will 
enhance the role of fluctuations which in turn are ex- 
pected to reduce the value of the transition temperature. 
Therefore, the qualitative properties of our paramagnetic 
solution may be relevant for real quantum spin glasses in 
their disordered phase. In this context two predictions 
that emerge from our work may provide useful insight 
for the analysis of experiments on quantum frustrated 
systems provided that T g is small enough, i) Mea- 
surements of the magnetic susceptibility as the transition 
is approached from above may indicate an anomalously 
low value of the effective Curie constant, ii) For T > T g 
a fraction of the spectral weight may be spread over a 
wide energy range and be difficult to distinguish from 
background noise. This fact, combined with the presence 
of a very strong and narrow central peak, may result in 
an apparent loss of spectral weight in neutron scattering 
experiments. 

Many interesting questions remain to be addressed, for 
instance, the effect of coupling the spin-system to an elec- 
tronic band. For a bandwidth larger than J one may 
expect that T g — > leading to the interesting physics of 
systems near quantum critical points uM. 
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